#### Online Appendix Material ####

library(tidyverse)
library(ggpubr)
library(patchwork)
library(synthdid)
library(ggrepel)
library(xtable)
library(kableExtra)
library(knitr)
library(stargazer)
library(modelsummary)
library(purrr)
library(gsynth)
library(plotrix)
library(sandwich)
library(lmtest)

set.seed(7390)

ward_full <- read.csv('data/ward_full.csv')

#### Creating Figure A.1 ####

hetro_plt <- ggplot(ward_full, aes(x = trnt))+
  geom_histogram(aes(y=..density..), position="identity", 
                 alpha=0.5, fill = 'red')+
  geom_density(colour = 'blue', alpha = .5)+
  facet_wrap(year ~ .)+
  scale_x_continuous('Turnout', labels = scales::percent, 
                     breaks = seq(0,1, by = 0.2))+
  scale_y_continuous('Density')+
  theme_minimal()+
  theme(panel.grid = element_blank(), 
        axis.text = element_text(size = 14),
        axis.title = element_text(size = 14),
        strip.text = element_text(size = 14))


ggsave('plots/hetro_plt.jpg', hetro_plt, 
       width = 35, height = 20, units = 'cm')

#### Creating Table A.1 ####

full_com_agg <- read.csv('data/full_com_agg.csv')

### Demo Props ###

full_com_agg <-  full_com_agg %>% 
  filter(loc != 'NA', 
         loc != 'City of London') %>%
  group_by(loc) %>% 
  arrange(loc) %>% 
  rename(Borough = 'loc',
         '\\% Black' = per_black,
         '\\% DE' = per_de,
         '\\% Retired' = per_ret, 
         '\\% Unemployed' = per_unemp, 
         '\\% Long-term Sick or Disabled' = per_dis)


demo_props <- knitr::kable(full_com_agg,
                           'pipe', digits = 2,
                           caption = '\\label{tab:demo_props} Proportions of each socio-demographic group by borough')

save_kable(demo_props, file = 'plots/demo_props.Rmd')

#### Creating Figure A.2 ####

ward_full <- ward_full %>% 
  filter(count == 5 & loc != 'Tower Hamlets') %>% 
  mutate(nbrom = ifelse(loc == 'Bromley', 1 , 0))

setup <- panel.matrices(as.data.frame(ward_full),
                        unit = which(names(ward_full) == 'ward_code'),
                        time =  which(names(ward_full) == 'year'),
                        outcome = which(names(ward_full) == 'trnt'),
                        treatment = which(names(ward_full) == 'treat'))

tau.hat <- synthdid_estimate(setup$Y, setup$N0, setup$T0)
tau.hat

tau.sc <- sc_estimate(setup$Y, setup$N0, setup$T0)
tau.sc

tau.did <- did_estimate(setup$Y, setup$N0, setup$T0)
tau.did

estimates <- list(tau.did, tau.hat, tau.sc)

names(estimates) <- c('Diff-in-Diff', 'Synthetic Diff-in-Diff', 'Synthetic Control')

print(unlist(estimates))

esti_plot <- synthdid_plot(estimates, facet.vertical=F,
                           control.name='Control', treated.name='Bromley',
                           lambda.comparable=TRUE, se.method = 'none',
                           trajectory.linetype = 1, line.width=2, 
                           effect.curvature=-.4, trajectory.alpha=.7, 
                           effect.alpha=1, diagram.alpha=1, 
                           onset.alpha=1) +
  theme(legend.position= 'bottom', 
        legend.direction='horizontal',
        legend.key=element_blank(), 
        legend.background=element_blank(),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14),
        panel.grid = element_blank(),
        axis.title = element_text(size = 14),
        axis.text = element_text(size = 14),
        strip.text = element_text(size = 14))+
  scale_color_discrete('Borough')+
  scale_y_continuous(labels = scales::percent, 
                     limits = c(0.25,0.75), expand = c(0,0),
                     breaks = seq(0.3, 0.7, by = 0.1))+
  scale_x_continuous(breaks = seq(2002, 2018, by = 4))+
  coord_cartesian(xlim = c(2002,2018))+
  xlab('Year')+
  ylab('Turnout')

esti_plot

ggsave('plots/esti_plot.png', esti_plot, 
       width = 35, height = 25, units = 'cm')

#### Creating Figure A.3 ####

vc.did.pl <- vcov(tau.did, method='placebo')
se.did.pl <- sqrt(diag(vc.did.pl))

vc.hat.pl <- vcov(tau.hat, method='placebo')
se.hat.pl <-  sqrt(diag(vc.hat.pl))

vc.sc.pl <- vcov(tau.sc, method='placebo')
se.sc.pl <- sqrt(diag(vc.sc.pl))

se.did.ja <-  sqrt(vcov(tau.did, method='jackknife'))
se.hat.ja <-  sqrt(vcov(tau.hat, method='jackknife'))
se.sc.ja <- sqrt(vcov(tau.sc, method='jackknife'))

se.did.boot <- sqrt(vcov(tau.did, 
                           method='bootstrap', 
                           replications = 1000))

se.hat.boot <-  sqrt(vcov(tau.hat, 
                            method='bootstrap', 
                            replications = 1000))

se.sc.boot <- sqrt(vcov(tau.sc, method='bootstrap', 
                          replications = 1000))

method <- c('Diff-in-Diff', 
            'Synthetic Difference-in-Difference', 
            'Syntehtic Control',
            'Diff-in-Diff', 
            'Synthetic Difference-in-Difference', 
            'Syntehtic Control',
            'Diff-in-Diff', 
            'Synthetic Difference-in-Difference', 
            'Syntehtic Control')

error_typ <- c('Bootstrap', 'Bootstrap', 'Bootstrap',
               'Placebo', 'Placebo', 'Placebo',
               'Jackknife', 'Jackknife', 'Jackknife')

minimum <- c(tau.did - 1.96*se.did.boot,
             tau.hat - 1.96*se.hat.boot,
             tau.sc - 1.96*se.sc.boot,
             tau.did - 1.96*se.did.pl,
             tau.hat - 1.96*se.hat.pl,
             tau.sc - 1.96*se.sc.pl,
             tau.did - 1.96*se.did.ja,
             tau.hat - 1.96*se.hat.ja,
             tau.sc - 1.96*se.sc.ja)

maxi <- c(tau.did + 1.96*se.did.boot,
          tau.hat + 1.96*se.hat.boot,
          tau.sc + 1.96*se.sc.boot,
          tau.did + 1.96*se.did.pl,
          tau.hat + 1.96*se.hat.pl,
          tau.sc + 1.96*se.sc.pl, 
          tau.did + 1.96*se.did.ja,
          tau.hat + 1.96*se.hat.ja,
          tau.sc + 1.96*se.sc.ja)

point <- c(tau.did, tau.hat, tau.sc,
           tau.did, tau.hat, tau.sc,
           tau.did, tau.hat, tau.sc)

error_dt <- data.frame(method, error_typ, minimum, maxi, point)

pla.boot.plt <- ggplot(error_dt, group = error_typ)+
  geom_errorbar(aes(ymin =  minimum, ymax = maxi,
                    x = method, colour = error_typ),
                alpha = 0.8, width = 0.25, size = 1.5,
                position = position_dodge(width = 0.5))+
  geom_point(aes(y = point,
                 x = method, colour = error_typ), 
             alpha = 0.8, size=2,
             position = position_dodge(width = 0.5))+
  coord_cartesian(ylim = c(-0.12, 0.04))+
  scale_y_continuous(breaks = seq(-0.12, 0.06, by = 0.02),
                     labels = seq(-0.12, 0.06, by = 0.02), 
                     'Treatment Effect')+
  scale_x_discrete('Method')+
  scale_colour_discrete('Error Estimation Method')+
  theme_minimal()+
  theme(panel.grid = element_blank(), 
        axis.text = element_text(size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14),
        axis.title = element_text(size = 14))+
  geom_hline(yintercept = 0, colour = 'red', linetype = 2)

pla.boot.plt 

ggsave('plots/pla.boot.plt.png', pla.boot.plt, 
       width = 35, height = 20, units = 'cm')

#### Creating Figure A.4 ####

ward_n2010 <- ward_full %>% 
  filter(year != 2010)

setup_n10 <- panel.matrices(as.data.frame(ward_n2010 ),
                            unit = which(names(ward_n2010) == 'ward_code'),
                            time =  which(names(ward_n2010) == 'year'),
                            outcome = which(names(ward_n2010 ) == 'trnt'),
                            treatment = which(names(ward_n2010) == 'treat'))

tau.hat_n10 <- synthdid_estimate(setup_n10$Y, setup_n10$N0,
                                 setup_n10$T0)
tau.hat_n10

#### Remove 2014 ####

ward_n2014 <- ward_full %>% 
  filter(year != 2014)

setup_n14 <- panel.matrices(as.data.frame(ward_n2014),
                            unit = which(names(ward_n2014) == 'ward_code'),
                            time =  which(names(ward_n2014) == 'year'),
                            outcome = which(names(ward_n2014) == 'trnt'),
                            treatment = which(names(ward_n2014) == 'treat'))

tau.hat_n14 <- synthdid_estimate(setup_n14$Y, setup_n14$N0,
                                 setup_n14$T0)
tau.hat_n14


#### Remove 2010 and 2014 ####

ward_n1014 <- ward_full %>% 
  filter(year != 2010 & year != 2014)

setup_n1014 <- panel.matrices(as.data.frame(ward_n1014),
                              unit = which(names(ward_n1014) == 'ward_code'),
                              time =  which(names(ward_n1014) == 'year'),
                              outcome = which(names(ward_n1014) == 'trnt'),
                              treatment = which(names(ward_n1014) == 'treat'))

tau.hat_n1014 <- synthdid_estimate(setup_n1014$Y, 
                                   setup_n1014$N0,
                                   setup_n1014$T0)
tau.hat_n1014




#### Tory Councils Only ####

ward_full_tories <- ward_full %>% 
  filter(loc == 'Barnet' | loc == 'Bromley' 
         | loc == 'Kensington and Chelsea' 
         | loc == 'Wandsworth' | loc == 'Westminster')

setup_tory <- panel.matrices(as.data.frame(ward_full_tories),
                             unit = which(names(ward_full_tories) == 'ward_code'),
                             time =  which(names(ward_full_tories) == 'year'),
                             outcome = which(names(ward_full_tories) == 'trnt'),
                             treatment = which(names(ward_full_tories) == 'treat'))

tau.hat_tory <- synthdid_estimate(setup_tory$Y, 
                                  setup_tory$N0,
                                  setup_tory$T0)
tau.hat_tory

### Error ###

vc.hat.pl_n10 <- vcov(tau.hat_n10, method='placebo')
se.hat.pl_n10 <- sqrt(diag(vc.hat.pl_n10))
se.hat.pl_n10 

vc.hat.pl_n14 <- vcov(tau.hat_n14, method='placebo')
se.hat.pl_n14 <- sqrt(diag(vc.hat.pl_n14))
se.hat.pl_n14

vc.hat.pl_n1014 <- vcov(tau.hat_n1014, method='placebo')
se.hat.pl_n1014 <- sqrt(diag(vc.hat.pl_n1014))
se.hat.pl_n1014

vc.hat.pl_tory <- vcov(tau.hat_tory, method='placebo')
se.hat.pl_tory <- sqrt(diag(vc.hat.pl_tory))
se.hat.pl_tory

pla.error_rm <- ggplot()+
  geom_errorbar(aes(ymin =  
                      c(tau.hat_n10 - 1.96*se.hat.pl_n10,
                        tau.hat_n14 - 1.96*se.hat.pl_n14,
                        tau.hat_n1014 - 1.96*se.hat.pl_n1014,
                        tau.hat_tory - 1.96*se.hat.pl_tory),
                    ymax = c(tau.hat_n10 + 1.96*se.hat.pl_n10,
                             tau.hat_n14 + 1.96*se.hat.pl_n14,
                             tau.hat_n1014 + 1.96*se.hat.pl_n1014,
                             tau.hat_tory + 1.96*se.hat.pl_tory),
                    x = c('2010 General Election Removed', 
                          '2014 EU Election Removed', 
                          'Both 2010 and 2014 Elections Removed',
                          'Only Tory Held Councils')), 
                colour = '#F8766D',
                alpha = 0.8, width = 0.25, size = 1.5)+
  geom_point(aes(y = c(tau.hat_n10, tau.hat_n14, 
                       tau.hat_n1014, tau.hat_tory),
                 x = c('2010 General Election Removed', 
                       '2014 EU Election Removed', 
                       'Both 2010 and 2014 Elections Removed',
                       'Only Tory Held Councils')), 
             colour = '#00BFC4', alpha = 1, size = 2)+
  coord_cartesian(ylim = c(-0.12, 0))+
  scale_y_continuous(breaks = seq(-0.12, 0.01, by = 0.03),
                     labels = seq(-0.12, 0.01, by = 0.03), 
                     'Treatment Effect')+
  scale_x_discrete('Design')+
  theme_minimal()+
  theme(panel.grid = element_blank(), 
        axis.text = element_text(size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14),
        axis.title = element_text(size = 14))+
  geom_hline(yintercept = 0, colour = '#F8766D', linetype = 2)

pla.error_rm 

ggsave('plots/pla.error_rm.png', pla.error_rm, 
       width = 35, height = 20, units = 'cm')
